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Abstract 

We study the random link travehng salesman problem, where lengths lij between 
city i and city j are taken to be independent, identically distributed random vari- 
ables. We discuss a theoretical approach, the cavity method, that has been proposed 
for finding the optimum tour length over this random ensemble, given the assumption 
of replica symmetry. Using finite size scaling and a renormalized model, we test the 
cavity predictions against the results of simulations, and find excellent agreement over 
a range of distributions. We thus provide numerical evidence that the replica sym- 
metric solution to this problem is the correct one. Finally, we note a surprising result 
concerning the distribution of fcth-nearest neighbor links in optimal tours, and invite 
a theoretical understanding of this phenomenon. 
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1 Introduction 



Over the past 15 years, the study of the travehng salesman problem (TSP) from the point of 
view of statistical physics has been gaining added currency, as theoreticians have improved 
their understanding of the relation between combinatorial optimization and disordered sys- 
tems. The TSP may be stated as follows: given sites (or "cities"), find the total length L 
of the shortest closed path ( "tour" ) passing through all cities exactly once. In the stochastic 
TSP, the matrix of distances separating pairs of cities is drawn randomly from an ensemble. 
The ensemble that has received the most attention in the physics community is the random 
link case, where the individual lengths Uj between city i and city j {i < j) are taken to be 
independent random variables, all identically distributed according to some p{l). The idea 
of looking at this random link ensemble, rather than the more traditional "random point" 
ensemble where cities are distributed uniformly in Euclidean space, originated with an at- 
tempt by Kirkpatrick and Toulouse to find a version of the TSP analogous to the earlier 
Sherrington-Kirkpatrick (SK) model for spin glasses. 

The great advantage of working with the random link TSP, rather than the (random 
point) Euclidean TSP, is that one may realistically hope for an analytical solution. A major 
breakthrough occurred with the idea, first formulated by Mezard and Parisi and later 
developed by Krauth and Mezard 0, that the random link TSP could be solved using the 
cavity method, an approach inspired by work on spin glasses. This method is based on 
assumptions pertaining to properties of the system under certain limiting conditions. The 
most important of these assumptions is replica symmetry. Although in the case of spin 
glasses, replica symmetry is violated 0, for the TSP there are various grounds for at least 
suspecting that replica symmetry holds The cavity solution then leads to a system of 

integral equations that can be solved — numerically at least — to give a prediction of the 
optimum tour length L in the many-city limit A^ ^ oo. 

In a previous article ||^, we have taken the random link distribution p{l) to match that of 
the distribution of individual city-to-city distances in the Euclidean case, using the random 
link TSP as a random link approximation to the Euclidean TSP. The approximation may 
seem crude since it neglects all correlations between Euclidean distances, such as the triangle 
inequality. Nevertheless, it gives remarkably good results. In particular, a numerical solution 
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of the random link cavity equations predicts large N optimum tour lengths that are within 
2% of the (simulated) (i-dimensional Euclidean values, for d = 2 and d = 3. In the limit 
d ^ oo, this gap shows all signs of disappearing. The random link problem, and its cavity 
prediction, is thus more closely related to the Euclidean problem than one might expect. 

The random link TSP is also, however, interesting in itself. Little numerical work has 
accompanied the analytical progress made — a shortcoming made all the more troubling by 
the uncertainties surrounding the cavity method's assumptions. In this paper we attempt to 
redress the imbalance, providing a numerical study of the finite size scaling of the random link 
optimum tour length, and arguments suggesting that the cavity solution is in fact correct. 
In the process, our numerics reveal some remarkable properties concerning the frequencies 
with which cities are connected to their fcth-nearest neighbor in optimal tours; we invite a 
theoretical explanation of these properties. 

2 Background and the cavity method 

In an attempt to apply tools from statistical mechanics to optimization problems, Kirk- 
patrick and Toulouse |^ introduced a particularly simple case of the random link TSP. The 
distribution of link lengths Itj was taken to be uniform, so that p{l) is constant over a fixed 
interval. In light of the random link approximation, one may think of this as corresponding, 
at large N, to the 1-D Euclidean case. (When cities are randomly and uniformly distributed 
on a line segment, the distribution of lengths between pairs of cities is uniform.) Although 
the 1-D Euclidean case is trivial — particularly if we adopt periodic boundary conditions, 
in which case the optimum tour length is simply the length of the line segment — the 
corresponding random link problem is far from trivial. 

The simulations performed by Kirkpatrick and Toulouse suggested a random link op- 
timum tour length value of Lju ^ 1.045 in the N ^ oo limit.Q Mezard and Parisi 
attempted to improve both upon this estimate and upon the theory by using replica tech- 
niques often employed in spin glass problems (for a discussion of the replica method in this 

"'^Here we work in units where the hne segment is taken to have unit length, and in order to match 
the normahzed 1-D Euchdean distribution, we let p{l) = 2 on [0,1/2). The 1-D Euclidean value, for 
comparison, would thus be — 1. Kirkpatrick and Toulouse, among others, choose instead p{l) = 1 on 
[0, 1), contributing an additional factor of 2 in Lijl which we omit when quoting their results. 
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context, see |^). This approach allowed them to obtain, via a saddle point approximation, 
many orders of the high-temperature expansion for the internal energy. They then extrap- 
olated down to zero temperature — corresponding to the global TSP optimum — finding 
L^L = 1.04 ± 0.015. Their analysis, like that of Kirkpatrick and Toulouse, was carried out 
only for the case of p{l) equal to a constant. 

Given the difficulties of pushing the replica method further, Mezard and Parisi then tried 
a different but related approach known as the cavity method [0]. This uses a mean-field 
approximation which, in the case of spin glasses, gives the same result as the replica method 
in the thermodynamic limit (A^ —>■ oo). As much of the literature on the cavity method 
has been prohibitively technical to non-specialists, we shall review the approach in more 
conventional language here, indicating what is involved in the case of the TSP. 

Both the replica and the cavity methods involve a representation of the partition function 
originally developed in the context of polymer theory |^, The approach consists of 

mapping the TSP onto an m-component spin system, writing down the partition function 
at temperature T, and then taking the limit m — 0. More explicitly, consider spins Sj, 
i = 1,...,N (corresponding to the A^ cities), where each spin Sj has m components Sf, 
a = 1, . . . , m, and where (Sj)^ = m for all i. The partition function is defined, in terms of a 
parameter uj, as 




L k<l 

where the integral is taken over all possible values spin values (the area measure is normalized 
so that / dfi{Sq) = 1), and Rij is related to the length lij between city i and city j as 
Rij = e^^^'''*''^/-^. Now employ a classic diagrammatic argument: let each spin product 
(Sa • Sb) appearing in the series be represented by an edge in a graph whose vertices are the 
A^ cities. The first-order terms (u) will consist of one-edge diagrams, the second-order terms 
(cu^) will consist of two-edge diagrams, and so on. What then happens when we integrate over 
all spin configurations? If there is a spin Sa that occurs only once in a given diagram, i.e., it 
is an endpoint, the spherical symmetry of Sa will cause the whole expression to vanish. The 
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non- vanishing summation terms in therefore correspond only to "closed" diagrams, where 
there is at least one loop. It may furthermore be shown that in performing the integration, 
any one of these closed diagrams will contribute a factor m for every loop present in the 
diagram If we then consider {Z — l)/m and take the limit m — 0, it is clear that only 
diagrams with a single loop will remain. Furthermore, since any closed diagram with more 
than links must necessarily contain more than one loop, only diagrams up to order 
will remain. Finally, take the limit uj (yo. The term that will then dominate in (|^) is the 
order term which, being a single loop diagram, represents precisely a closed tour passing 
through all cities. We may write it without the combinatorial factor A^! by expressing it 
as a sum over ordered pairs in the tour, and we thus find: 

Z — 1 

[J— >oo TV-link single loops 

(«1,«2, 

Af-city tours 

where L is the total tour length. What we obtain is exactly the partition function for the 
traveling salesman problem, with the correct canonical ensemble Boltzmann weights, using 
the tour length as the energy to be minimized (up to a factor N^/'^^ necessary for the energy 
to be extensive). 

The idea behind the cavity method is as follows. Since all spin couplings Rij in (|1|) 
are positive (ferromagnetic), we expect the m-component spin system to have a non-zero 
spontaneous magnetization in equilibrium. Now add an (A^ + l)th spin to the system; it too 
acquires a spontaneous magnetization (Sat+i). Let us obtain the thermodynamic observables 
of the new system (in particular (Sat+i) itself) in terms of the earlier magnetizations (Sj)' 
from before the (A^ + l)th spin was added — hence the notion of a "cavity". 

In order to compute these relations, an important mean-field assumption is made: that at 
large A^, any effect spin A^ + 1 feels from correlations among the A^ other spins is negligible. 
We justify this in the following way. Although all spins in (|l|) are indeed coupled, the coupling 
constants Rij decrease exponentially with length Zj^, and so effective interactions arise only 
between very near neighbors. But a crucial property of the random link model is that the 
near neighbors of spin A^ + 1 are not generally near neighbors of one another: they are near 
neighbors of one another only with probability 0(1/A^). Thus, when considering quantities 
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involving spin + 1, the effect of direct interactions between any two of its neiglibors is 
0{1/N), and decays to zero in tlie limit N ^ oo. We tlierefore replace (|I|) by the mean-field 
partition function 

. N . N N+1 

Zmf = / n djJ-i^q) / dn{SN+i) exp{uj J^^i^N+i Si ■ Sat+i + ^ Si ■ hi). (5) 

q=l 1=1 i=l 

By definition, if spin + 1 were removed, we would recover the "cavity magnetizations" 
(Si)'. This requirement is sufficient to specify the fields hj. Stripping out spin + 1 from 
(|^) leaves us simply with a product of integrals Yiql dfi{Sg) exp(Sg ■ hg), whose logarithmic 
derivative with respect to hj must then give the magnetization (Sj)'. We may obtain this 
expression by expanding the integrands, taking advantage of the identity / d^{Si) S^S^ = dap 
for all spin components a and (3, as well as the nilpotency property |^ that in the limit 
m — s> 0, integrating the product of more than two components of Sj gives zero. (This is 
analogous to the property used earlier in the diagrammatic expansion.) The result is 

<^'>' = l + (^)V2 - 

Note that this specifies hi for 1 < i < A^; hjv+i has been introduced purely for analytical 
convenience, and will ultimately be set to 0. 

Without loss of generality, let us assume the spontaneous magnetizations of the system 
to be directed exclusively along component 1. This may be imposed, for instance, by ap- 
plying an additional infinitesimal field directed along component 1. Physically, however, 
the assumption that distant spins are uncorrelated also means that this infinitesimal field 
is sufficient to select a single phase or equilibrium state, thus giving rise to a unique ther- 
modynamic limit. From the point of view of dynamics, a consequence is that two "copies" 
of the system will evolve to the same equilibrium distribution. This property is known as 
replica symmetry, and has been central to the modern understanding of disordered systems.^ 
Replica symmetry is in fact known to be broken in spin glasses; if one uses, for instance, the 
replica symmetric solution of the SK model, one finds a ground state energy prediction that 
is inaccurate by about 5% |p. However, this does not mean that replica symmetry breaking 
occurs in all related problems of high complexity (the TSP and the spin glass both fall into 
^An analogous property was used to obtain the replica solution mentioned earlier. 
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the NP-hard class of computational complexity). Showing that the (replica symmetric) cav- 
ity solution correctly predicts macroscopic quantities for the random link TSP would suggest 
that the TSP, unlike a spin glass, does indeed exhibit replica symmetry. 

In order to obtain the cavity solution, consider the mean- field expression (^. Taking 
advantage of nilpotency, as well as the fact that Rij is non-negligible only with probability 
0{1/N), we may expand @ and obtain in the large limit: 



X 



. (hAT+l) ^ CU-Rj^AT+i hj ■ Hat+i ^ LlJ Rj^N+lRk,N+l^j -^k 

~^ 9 ^ ^ 1 _L (U .\2 /O ^ 2^ 



1 + (h,)V2 i<,^<^ [1 + (h,)V2][l + (h.)V2]_ 

Differentiating Zjv/f with respect to hjv+i and then setting h^v+i = yields an expression 
for (Sat+i) in terms of the remaining hj, or equivalently in terms of the cavity magneti- 
zations (Sj)'. This expression simplifies further simplifies at large uj. Recalling that the 
magnetization is by construction directed along component 1, we obtain 0: 

^ J2l<j<k<N Rj,N+lRk,N+l{Sj)' (Si)' 

(The factor u in the denominator may be avoided, if need be, by applying a uniform rescaling 
factor y/uj to all magnetizations.) Thus, using the mean-field approach, we can express the 
magnetization of the (A^ + l)th spin in terms of what the other magnetizations would be in 
the absence of this (A^ + l)th spin. 

While these quantities have been derived for a spin system whose partition function is 
given by Z, we are interested in the TSP whose partition function is given by (^). Consider 
an important macroscopic quantity for the TSP: the frequency with which a tour occupies 
a given link. Define riij to be 1 if the link ij is in the tour, and otherwise. Since the total 
tour length (energy) is L = Y^i^jfiijlij, the mean occupation frequency (rijj), averaged over 
all tours with the Boltzmann factor, is simply found from the logarithmic derivative of (^) 
with respect to kj. Using Zmf — 1 in place of Z — 1, and proceeding as above, we obtain in 
the limit u ^ oo: 

i„ \_D T.j^iRj,N+i{S])' 

\ni,N+i) - Ki^N+i{^i ) ^ B B /Qi\/ /on/- 

2^1<j<k<N -'^j,N+l-n'k,N+l\^j / Wfc/ 
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The relations and have been derived for a single realization of the Rij^s. In 
the ensemble of instances we consider here, the thermal averages become random variables 
with a particular distribution. As far as (j^) is concerned, we may treat the magnetizations 
{S})' as independent identically distributed random variables. Furthermore, the existence 
of a thermodynamic limit in the model requires that at large A^, have the same 

distribution as the cavity magnetizations; this imposes, for a given link length distribution 
p{l), a unique self-consistent probability distribution of the magnetizations. From (P), one 
can then find the probability distribution of (n7v+i,i), and in turn, taking the T — > limit, 
the distribution V{1) of link lengths / used in the optimal tour (at A^ —>■ cxd). 

Krauth and Mezard Q carried out this calculation, for p(/) corresponding to that of the 
d-dimensional Euclidean case, namely 

Of course, pd{l) must be cut off at some finite / in order to be normalizable; precisely how 
this is done is unimportant, however, since only the behavior of pd(/) at small / is relevant 
for the optimal tour in the N ^ oo limit. The result of Krauth and Mezard's calculation is: 

where H^i^x) is the solution to the integral equation 

From Vd{l), one may obtain the mean link length in the tour, and thus the cavity prediction 
L^^ for the total length of the tour. Introducing the large A^ asymptotic quantity Pmid) = 
lim^^oo Lrl{N, d) /N^~^^'^, the cavity prediction P^^id) is then: 

P^^id) = lim N'/" / lVd{l)dl 

rj r + co 

= -/ Hd{x)[l + Hd{x)]e~''-^^Ux. (13) 



At d = 1, Krauth and Mezard solved these equations numerically, obtaining /^^^(l) = 
1.0208 ... It is difficult to compare this with Kirkpatrick's value of Prl{^) ~ 1.045 from direct 
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simulations (as no error estimate exists for the latter quantity), however an analysis |n] of 



recent numerical results by Johnson et al. [O gives /?/?/,(!) = 1.0209±0.0002, lending strong 



credence to the cavity value. Krauth and Mezard also performed a numerical study of Vi{l). 
They found the cavity predictions to be in good agreement with what they found in their 
own direct simulations. Further numerical evidence supporting the assumption of replica 
symmetry was found by Sourlas [§, in an investigation of the low temperature statistical 
mechanics of the system. Thus, for the 1^ distribution at d = 1, there is good reason to 
believe that the cavity assumptions are valid and that the resulting predictions 
large N, so that /^^^(l) = 



At higher dimensions, the values of (3'^i{d) were given by the present authors in |T3[, and 
a large d power series solution for /^^^((i) was derived 0]: 



\ 2 - In 2 - 27 „ 



(14) 



where 7 represents Euler's constant (7 = 0.57722. . .). But is the cavity method exact — 
that is, is /3^i((i) = PRL^d) — for all d, or is d = 1 simply a pathological case (as it 
is in the Euclidean model, where /?£;(!) = 1)? While it appears sensible to argue that the 
qualitative properties of the random link TSP are insensitive to d, there is as yet no evidence 
that replica symmetry holds for d ^ 1. Our purpose here is to provide such evidence by 
numerical simulation, as has been done, for instance, in a related combinatorial optimization 
problem known as the matching problem KM. We now turn to this task, considering 



first the d = 2 case, and then a "renormalized" random link model that enables us to verify 



numerically the 0{l/d) coefficient predicted in (p!^). 



3 Numerical analysis: d = 2 case 

We have implicitly been making the assumption so far, via our notation, that as — oo 
the random variable LjiL(N,d)/N^~^^'^ approaches a unique value PuLid) with probability 
1. This is a property known as self-averaging. The analogous property has been shown for 



the Euclidean TSP at all dimensions |16|. For the random link TSP, however, the only case 



where a proof of self- averaging is known is in the d —>■ oo limit, where a converging upper 
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Figure 1: Distribution of 2-D random link rescaled tour length — {Lrl)) / for increasing 

values of N. Plus signs show = 12 (100,000 instances used), squares show = 17 (100,000 

instances used), diamonds show = 30 (4,000 instances used), and dots show = 100 (1,200 
instances used). Solid lines represent Gaussian fits for each value of A^ plotted. 



and lower bound give in fact the exact result |1T7| : 



PRLid) 



d 



27re 



1 + 



(15) 



Comparing this with ([T^), we may already see that P^^id) ~ PRL^d) when d ^ oo, and so 
the cavity prediction is correct in the infinite dimensional limit. 

For finite d, however, it has not been shown analytically that (3RL{d) even exists. To 
some extent, the difficulty in proving this can be traced to the non-satisfaction of the trian- 
gle inequality. The reader acquainted with the self-averaging proof for the Euclidean TSP 
may see that the ideas used there are not applicable to the random link case; for instance, 
combining good subtours using simple insertions will not lead to near-optimal global tours, 
making the problem particularly challenging. Let us therefore examine the distribution of 
d = 2 optimum tour lengths using numerical simulations, in order to give empirical support 
for the assertion that the N —>■ oo limit is well-defined. 
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The algorithmic procedures we use for simulations are identical to those we have used 
in an earlier study concerning the Euclidean TSP for details, the interested reader is 
referred to that article. Briefly stated, our optimization procedure involves using the LK 
and CLO local search heuristic algorithms [IE, rU| where for each instance of the ensemble we 
run the heuristic over multiple random starts. LK is used for smaller values of (A^ < 17) 
and CLO, a more sophisticated method combining LK optimization with random jumps, for 
larger values of (A^ = 30 and A^ = 100). There is, of course, a certain probability that 
even over the course of multiple random starts, our heuristics will not find the true optimum 
of an instance. We estimate the associated systematic bias using a number of test instances, 
and adjust the number of random starts to keep this bias at least an order of magnitude 
below other sources of error discussed below. (At its maximum — occurring in the A^ = 100 
case — the systematic bias is estimated as under 1 part in 20,000.) 

Following this numerical method, we see from our simulations (Figure |I]) that the distri- 
bution of LfiL{N, 2)/\/N becomes increasingly sharply peaked for increasing A^, so that the 
ratio approaches a well-defined limit (3rl{2). Furthermore, the variance of Ljil{N, 2) remains 
relatively constant in A^ (see Table |I|), indicating that the width a for the distribution shown 
in the figure decreases as strongly suggesting a Gaussian distribution. Similar results 

were found in our Euclidean study (albeit in that case with a being approximately half of its 
random link value). This is precisely the sort of behavior one would expect were the central 
limit theorem to be applicable. 

Let us now consider the large A^ limit of Lrl{N, 2)/^/N, as given by numerical simula- 
tions. In the Euclidean case, it has been observed M that the finite size scaling law can be 



Table 1: Variance of the non-rescaled optimum tour length Lrl{N, 2) with increasing A^. 

N (7^ ^ instances used 

~l2 0.3200 100,000 

17 0.3578 100,000 

30 0.3492 4,000 

100 0.3490 1,200 
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Figure 2: Finite size scaling of mean optimum tour length for d = 2. Best fit (x^ = 4.46) is given 
by: (Li^z,(iV,2))/iVV2 = 0.7243(1 + 0.0322/Af - l.SSe/iV^)^ Error bars show one standard deviation 

(statistical error). 



written in terms of a power series in 1/A^. The same arguments given there apply to the 
random link case, and so we may expect the ensemble average {L^l^N ,2)) to satisfy 

{LRL{N,d))=(3RL{d)N^-^''' 

In order to obtain {Lb,l{N ^2)) at a finite value of from simulations, we average over a 
large number of instances to reduce the statistical error arising from instance-to-instance 
fluctuations. Figure |^ shows the results of this, with accompanying error bars, fitted to the 
expected finite size scaling law (truncated after 0{1/N'^)). The fit is a good one: = 4.46 for 
5 degrees of freedom. As in 0, we may obtain an error estimate on /5rl(2) by noting that if 
we take the extrapolated value and add or subtract one standard deviation, and then redo the 
fit with this as a fixed constant, will increase by 1. We thus find Prl{2) = 0.7243±0.0004, 
in very good agreement with the cavity result of (3'^j^{2) = 0.7251 . . . The discrepancy between 
the two is consistent with the statistical error (two standard deviations apart), and in relative 



1 + 



A{d) 
N 



+ 



(16) 
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terms is approximately 0.1%. The fit in Figure furthermore, appears robust with respect 
to sub-samples of the data; even if we disallow the use of the = 100 data point in the fit, 
the resulting asymptotic value is still within 0.25% of the cavity prediction. By comparison, 
recall that the error in the replica symmetric solution to the SK spin glass ground state 
energy is of the order of 5% f^. 

Another quantity that Krauth and Mezard studied in their d = 1 numerical investigation 
Q was the optimum tour link length distribution Vd{l) given in (pH]) . Let us consider 
V2{1), and following their example, let us look specifically at the integrated distribution 
/^(Z) = J^Vdil') dl' . The cavity result for Idil) can, like (^mid), be computed numerically to 
arbitrary precision. In Figure ^ we compare this with the results of direct simulations, for 
d = 2, at increasing values of A^. The improving agreement for increasing A^ (within 2% at 
A^ = 100) strongly suggests that the cavity solution gives the exact A^ — oo result. 

Finally, it is of interest to consider one further quantity in the d = 2 random link simu- 

1.0 
0.8 
0.6 
0.4 
0.2 
0.0 

0.0 0.5 1.0 1.5 2.0 

/ 

Figure 3: Integrated probability distribution of link lengths in the optimal tour, for d = 2, using 
rescaled length I = ly/N . Plus signs represent = 12 simulation results, dots represent A^ = 100 
simulation results, and solid line represents cavity prediction. 
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Figure 4: Frequencies with which /cth-nearest neighbors are used in optimal 2-D random link tours. 
Plus signs show values for N = 12, squares for = 17, diamonds for A'^ = 30, and dots for N = 100. 
Best exponential fit (straight line on log plot) is shown for = 100 data. 

lations, for which there is at present no corresponding cavity prediction: the frequencies of 
"neighborhood rank" used in the optimal tour, that is, the proportion of hnks connecting 
nearest neighbors, 2nd-nearest neighbors, etc. Sourlas P] has noted that in practice in the 
d = 1 case, this frequency falls off rapidly with increasing neighborhood rank — suggest- 
ing that optimization heuristics could be improved by preferentially choosing links between 
very near neighbors. Our simulations show (see Figure ^) that for d = 2 the decrease is 
astonishingly close to exponential. We may offer the following qualitative explanation for 
this behavior. An optimal tour will always try to use links to the closest neighbors possible. 
While the constraint of a closed loop may force it in rare cases to use neighbors of high rank, 
this will apply only to a very small number of links in the tour. Connecting a point to, say, 
its fcth-nearest neighbor will for the most part be profitable only when this neighbor is not 
much further away than the k — 1 nearer neighbors. In other words, the lengths from the 
point to its k closest neighbors would have to be nearly degenerate. Since a fc-fold degen- 
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eracy of this sort is the product of — 1 unhkely events, it is in fact quite natural that the 
probabihty of such an occurrence is exponentially small in k. 

We therefore conjecture that the neighborhood frequency function will fall exponentially 
in k at large k. We expect this behavior to hold in any dimension, and for that matter, in 
the Euclidean TSP as well. Similar and even stronger numerical results have been reported 
pT| in another link-based combinatorial optimization problem, the matching problem. An 
analytical calculation of the neighborhood frequency may indeed turn out to be feasible using 
the cavity approach, thus providing a theoretical prediction to accompany our conjecture. 
We consider this a significant open question. 

4 Numerical analysis: renormalized model 

In this section we will consider a different sort of random link TSP, proposed in allowing 
us to test numerically the 1/d coefficient predicted by the cavity result (Jji^). The approach 
involves introducing a mapping that shifts and rescales all the lengths between cities. By 
taking the limit d oo, one obtains a rf- independent random link model having an expo- 
nential distribution for its link lengths. This "renormalized" model was outlined in [0; we 
present it here in further detail. We then perform a numerical study of the model, which 
enables us to determine the large d behavior of the standard d- dimensional random link 
model. 

Let us define {Di{N,d)) to be the distance between a city and its nearest neighbor, 
averaged over all cities in the instance and over all instances in the ensemble.0 For large d, 
it may be shown that 



where 7 is Euler's constant. It is not surprising that this quantity is reminiscent of ([T5|) , 
since N^/^{Di(N, d)) represents precisely a lower bound on Pniid). 

In order to obtain the renormalized model, consider a link length transformation mak- 
ing use of {Di{N,d)). For any instance with link lengths kj (taken to have the usual 

•^Note that {Di[N,d)) itself does not involve the notion of optimal tours, or tours of any sort for that 
matter. 




(17) 
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distribution ( ]TD| ) corresponding to d dimensions), define new link lengths Xij = d\lij — 
{Di{N,d))]/ {Di{N,d)). The Xij are "lengths" only in the loosest sense, as they can be 
both positive and negative. The optimal tour in the Xij model will, however, follow the same 
"path" as the optimal tour in the associated lij model, since the transformation is linear. Its 
length Lx{N, d) will simply be given in terms of LjulN, d) by: 



r^.r. ^ LnL{N,d)-N{D^{N,d)) 
Lx[i\,d) = d , „ , ,^ , , so 



LRL{N,d) = N{D,{N,d)) 



{DiiN,d)) 

L^{N,d) 



1 + 



dN 



(18) 
(19) 



In the standard dimensional random link model, PuLid) = limN^ao Lrl{N, d) /N^ so 

L,{N,dy 



PRLid) 



lim N^''^{Di{N,d)) 



1 + 



dN 



27ie 



(nd) 



l/2d 



7 



d 



d^ 



lim 

N^co 



and at large d, using (17), 

L,{N,dy 



1 + 



dN 



(20) 



As PRL^d) is a well-defined quantity, there must exist a value fi{d) such that 
lim^^oo L,{N, d)/N = fi{d). 

Now, what will be the distribution of "renormalized lengths" p{x) corresponding to p(0? 
Prom (p!0| ) and the definition of the Xy, 



p{x) 



T{d/2 + l) d 



-, and substituting for I, 



d/2 



d-1 



{D,{N,d)y 



T{d/2 + l) 

In the limit — > 00, we thus obtain from (0) the large d expression: 



(21) 



p{x) 



d/2 



T{d/2 + l] 



i + ^V"iv- 



d \ 



d/2 



dj 



27ie 



iV-Ml--j + . 



X 



d-1 



by Stirling's formula 



1 + 01- 



(22) 



In the limit d 00, p{x) will be independent of d; the same must then be true for Lx{N, d), 
and consequently for p{d). 
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Let us now define the renormalized model as being made up of link "lengths" Xij in 
this limit. This results in a somewhat peculiar random link TSP, no longer containing the 
parameter d. Its link length distribution is given by the d oo limit of (P^), 



p(x) = ^ exp(x — 7), 



(23) 



and its optimum tour length satisfies 



lim 



N 



/i, 



(24) 



where we have dropped the d argument from these (now c?- independent) quantities. By 
performing direct simulations using the distribution (^) — cut off beyond a threshold value 
of X, as was done for pd{l) — we may find the value of fi numerically. 

Finally, let us relate this renormalized model to the standard (i-dimensional random link 




-0.5 0.0 0.5 
Rescaled tour length 



1.5 



Figure 5: Distribution of renormalized random link rescaled tour length {L^ — {Lx))/N for increas- 
ing values of N. Plus signs show = 12 (100,000 instances used), squares show = 17 (100,000 
instances used), diamonds show = 30 (4,000 instances used), and dots show A^ = 100 (1,200 
instances used). Solid lines represent Gaussian fits for each value of A^ plotted. 
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0.722 



0.02 



0.04 0.06 



0.08 



0.1 



Figure 6: Finite size scaling of renormalized model optimum. Best fit {x^ = 5.23) is given by: 
{L^{N))/N = 0.7300(l+0.3575/iV-2.791/iV2)^ Ej-j-or gj^Q^ Q^e standard deviation (statistical 



error I. 



model. In light of (^4[), we may rewrite (pOD and obtain the result given in 0: 



I d 
27re 



{nd) 



1/2(1 



1 + 



d 



(25) 



The value of fi in the renormalized model therefore gives directly the 1/d coefficient for the 
(non-renormalized) PRiid). 

We now carry out these direct simulations for the renormalized model. Figures ^ and ^ 
show our numerical results. In Figure ^ we see that just as in the d = 2 case, the distribution 
of the optimum tour length becomes sharply peaked at large and the asymptotic limit fi 
is well-defined. Via (p5D , this provides very good reason for believing that (3RL{d) is well- 
defined for all d, and that self-averaging holds for the random link TSP in general. In Figure 
^ we show the finite size scaling of {Lx{N)/N) . The fit is again quite satisfactory (with 
X^=5.23 for 5 degrees of freedom), giving the asymptotic result /i = 0.7300 ± 0.0010. The 
resulting value for the 1/d coefficient in PRL^d) is then /i — 7 = 0.1528 ± 0.0010, in excellent 
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10° 




^ I , , , , I , , , , I , , , , I , , , , I , , , , I , , , , I , , , , I , , , , I , , , 

123456789 10 
Neighborhood rank (k) 

Figure 7: Frequencies with which fcth-nearest neighbors are used in optimal renormalized random 
link tours. Plus signs show values for A'^ = 12, squares for N = 17, diamonds for N = 30, and dots 
for = 100. Best exponential fit (straight line on log plot) is shown for N = 100 data. 

agreement (error under 0.3%) with the cavity prediction 2 — In 2 — 27 = 0.1524 . . . given in 

(HI)- 

Again, as in the d = 2 case, let us briefly consider the frequencies of /cth-nearest neighbors 
used in optimal tours. These frequencies are given in Figure |^ for the renormalized model. 
Even though the exponential flt is not as good as in the d = 2 case, it is still striking 
here. What does this tell us, in turn, about the standard random link TSP? Recall that the 
renormalized model arises from the d ^ 00 limit of the (i-dimensional (non-renormalized) 
model, and that the mapping ([I8D preserves the optimum tour for any given instance. These 
/cth-neighbor frequency results are thus the d ^ 00 limiting frequencies for the ci-dimensional 
random link TSP (and most likely for the Euclidean TSP also). This gives further support 
to our conjecture that the exponential law holds for all d, and suggests as a consequence 
that the "typical" neighborhood rank k used in optimal tours remains bounded for all d. 
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5 Conclusion 

The random link TSP has interested theoreticians primarily because of its analytical 
tractability, allowing presumably exact results that are not possible in the more traditional 
Euclidean TSP. Other than in the d = 1 case, however, it has attracted little attention. In 
this paper we have provided a numerical study of the random link TSP that was lacking up 
to this point, addressing important unanswered questions. Through simulations, we have 
tested the validity of the theoretical predictions derived using the cavity method. While in 
other disordered systems, such as spin glasses, the replica symmetric solution gives values 
of macroscopic quantities that are inexact (typically by several percent), in the random link 
TSP it shows all signs of being exact. We have studied various link-based quantities at d = 2 
and found that the numerical results confirm the cavity predictions to within 0.1%. Further- 
more, we have confirmed, by way of simulations on a renormalized random link model, that 
the analytical cavity solution gives a large d expansion for the optimum tour length whose 
1/d coefficient is correct to within well under 1%. The excellent agreement found at d = 1 
[^, 0, d = 2, and to 0{l/d) at large d, then suggest strongly that the cavity predictions are 
exact. This provides indirect evidence that the assumption of replica symmetry — on which 
the cavity calculation is based — is indeed justified for the TSP. 

Finally, our random link simulations have pointed to a surprising numerical result. If 
one considers the links in optimal tours as links between fcth-nearest neighbors, at = 2 the 
frequency with which the tour uses neighborhoods of rank k decreases with k as almost a 
perfect exponential. Encouraged by similar results in the renormalized model, we conjecture 
that this property holds true for all d, as well as in the Euclidean TSP. As no theoretical 
calculation presently explains the phenomenon, we would welcome further investigation along 
these lines. 
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